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ABSTRACT 


An  analytical  model  was  formulated  for  the  calculation  of 
unsymmetrical  pressure  distributions  inside  the  cavities  of  a 
multi-stage  labyrinth  seal.  The  model  is  based  on  a  set  of  linearized 
continuity  and  momentum  equations  for  the  throughflow,  and  can  account 
for  the  effects  of  rotation,  whirl,  flow  swirl,  differential  gap  width 
and  (limited)  compressibility.  Preliminary  calculations  show  excellent 
agreement  with  non-whirling  side-force  data  in  the  literature.  Also, 
models  are  being  constructed  that  will  permit  calculation  of  the  effect 
of  these  seal  forces  on  the  running  stability  of  flexible,  highly 
loaded  turbomachines,  with  the  Space  Shuttle  Main  Engine  serving  as  the 
prototype  of  such  machines. 
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FORWARD 


This  report  describes  work  done  at  the  Gas  Turbine  and  Plasma 
Dynamics  Laboratory  (GT&PDL)  at  the  Massachusetts  Institute  of 
Technology  for  the  Air  Force  under  Grant  No.  AFOSR-83-0034.  Dr. 
Anthony  K.  Amos  was  the  technical  monitor. 

The  work  reported  herein  was  performed  during  the  period  1 
December  i982  through  30  November  1983.  Two  graduate  students,  Otto 
W.K.  Lee  and  Ya-Pei  Chang,  worked  under  the  supervision  of  Professors 
Manuel  Martinez-Sanchez  and  John  Dugundji  during  this  time. 


dynamic-structural  interaction  of  labyrinth  seals.  This  research  is  in 
an  area  of  limiting  technology  in  modern  gas  turbine  engines.  A  better 
understanding  of  labyrinth  seals  and  their  improved  design  would  reduce 
stage  and  secondary  flow  leakages,  and  thus  lead  to  increased  engine 
efficiency  and  performance.  This  improvement  should  be  accomplished 
without  compromising  the  rotor  dynamic  stability,  which  is  known  to  be 
sensitive  to  the  fluid  dynamic  forces  arising  at  the  seals  (see 
References  (1)  and  (2)). 

The  present  investigation  was  divided  into  two  main  efforts,  one 
involving  a  fluid  dynamics  investigation  of  the  seal  forces  developed 
in  a  rotating  labyrinth  seal,  and  another  involving  the  rotordynamic 
and  structural  dynamic  stability  of  a  rotating  shaft  system  in  the 
presence  of  labyrinth  seals  and  conventional  bearings.  The 
investigation  was  envisaged  as  extending  over  more  than  one  year. 


WORK  PERFORMED 

(a)  Fluid  Dynamics  of  Seals 

The  literature  for  fluid  dynamic  modeling  of  labyrinth  seals  was 
reviewed.  Recent  efforts  at  comprehensive  modeling  of  labyrinth  seals 
have  been  made  by  Kostyuk' 3 1 • ( 4 1  and  Iwatsubo' 5 1 .  These  models 
considered  multicavity  seals  with  the  shaft  executing  orbiting  motion 
without  or  with  simultaneous  canting' 3 1 ' 1 4 1 .  Both  axial  and 


rotational  fluid  motion  was  considered  in  the  glands,  although  quite 


simplified  equations  were  used  to  describe  the  flow  in  each  cavity  and 
past  each  tooth.  Numerical  techniques  were  then  resorted  to  in  order 
to  obtain  overall  forces  from  the  original  partial  differential 
equations  in  time  and  azimuthal  angle. 

In  parallel  with  these  analytical  attacks,  detailed  pressure 
measurements  inside  oscillating  multicavity  seals  have  also  become 
recently  available  (Ref.  6).  These  data  address  a  wide  range  of 
variables  of  interest,  including  relative  rotational/whirl  speeds, 
amplitudes,  pre-swirl,  number  of  cavities  per  seal,  etc.  They 
generated  both  pressure  distributions  in  the  glands  and  integrated 
spring  and  damping  coefficients,  for  both  colinear  and  quadrature 
forces. 

The  methods  of  Refs  (3)  to  (5)  held  promise  of  providing  a 
comprehensive  and  practical  model  capable  of  quantitative  predictions 
at  the  design  stage.  Accordingly,  a  model  of  fluid  forces  in  a 
labyrinth  seal  was  fashioned  after  that  of  Iwatsubo' 5 1  .  Basically,  it 
consisted  of  a  momentum  and  a  continuity  equation  for  each  of  a  number 
of  seal  glands  executing  small  elliptic  satellite  oscillations.  The 
inter-cavity  flows  were  calculated  using  constant  discharge 
coefficients  and  accounting  for  compressibility  effects;  the  flow 
azimuthal  velocity  was  calculated  from  a  balance  between  friction  on 
the  fixed  and  movable  parts  of  the  seal,  using  conventional  turbulent 
friction  coefficients.  After  linearization  for  small  perturbations, 
there  resulted  a  set  of  linear,  periodic  coefficient  equations  for 
perturbation  pressures,  velocities,  etc.,  in  each  of  the  glands  of  the 


seal.  To  the  first  approximation,  harmonic  balance  was  used  to  reduce 
this  to  a  system  of  partial  differential  equations  in  time  and  azimuth, 
with  constant  coefficients.  The  direct  and  quadrature  components  of 
each  variable  were  then  conveniently  isolated,  and  by  restricting 
attention  to  the  quasi-steady  state  forced  state  (growth  or  decay  time 
much  longer  than  a  period) ,  the  solution  reduced  to  that  of  an 
algebraic  set  of  equations  in  the  respective  direct  and  quadrature 
amplitudes  and  phase  angles  of  all  the  variables.  Integration  around 
the  seal  then  related  these  amplitudes  to  the  desired  force 
coefficients. 

The  method  described  above  is  a  marked  improvement  over  those  used 
in  Ref.  (5)  for  example,  where  finite  difference  approximations  were 
used  to  solve  the  governing  equations.  It  is  anticipated  that  this 
will  permit  introduction  of  refinements  to  the  model  while  still 
remaining  computationally  economical. 

The  above  labyrinth  seal  formulation  is  being  written  up  as  a 
Master’s  thesis  by  one  of  the  two  graduate  students,  and  will  become 
available  in  February  1984,  Ref.  7.  Also,  a  paper  based  on  this 
formulation  will  be  presented  at  the  forthcoming  3rd  Workshop  on  Rotor 
Dynamic  Instability  Problems  in  High  Performance  Turbomachinery  (Texas 
A&M,  May  28-30,  1984).  A  copy  of  this  paper,  giving  the  mathematical 
formulation  of  the  labyrinth  seals  model,  is  included  here  as  Appendix 
A  of  this  report.  Also,  the  Master’s  Thesis  by  A.  Celorio  (Ref.  9)  was 
instrumental  in  initiating  our  seal  modeling  effort. 


(b)  Rotordynamic  and  Structural  Dynamic  Stability 

Some  basic  literature  for  the  dynamic  stability  of  rotor  systems 
was  reviewed  (Ref  10-14)  together  with  some  more  recent  work  related  to 
the  space  shuttle  fuel  turbopump  (Refs.  15,  16).  Using  the  traditional 
linear  methods  of  rotor  analysis,  a  simple  model  of  the  high  pressure 
fuel  turbopump  (HPFTP)  for  the  space  shuttle  is  being  set  up  based  on 
the  structural  dynamic  data  given  by  Muller  in  Ref.  16.  Equations  of 
motion  have  been  written  for  the  rotating  shaft  assembly  on  the 
flexible  supporting  case  structure,  including  provisions  for  labyrinth 
seal  forces.  The  labyrinth  seal  forces,  as  determined  from  the  fluid 
dynamic  part  of  the  investigation,  generally  involve  frequency 
dependent  stiffness  and  damping  terms.  These  seal  forces  are  being 
fitted  here  by  means  of  Pade  approximations  to  give  frequency 
independent  (constant  coefficient)  stiffness  damping,  mass  and  lag 
terms  to  be  used  in  a  traditional  eigenvalue  analysis  to  examine 
possible  rotor  instabilities.  This  is  similar  to  procedures  used  in 
representing  aerodynamic  forces  in  aeroelastic  analysis  of  aircraft 
(Ref.  17). 

The  above  rotordynamic  formulation  is  currently  being  pursued  by 
the  other  of  the  two  graduate  students,  and  will  become  available  as  a 
Master’s  Thesis  during  the  summer  of  1984. 
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APPENDIX  A 


Copy  of  paper  to  be  presented  at  the  3rd  Workshop  on  Rotordynamic 
Instability  Problems  in  High  Performance  Turbomachinery,  Texas  A&M 


Univ.  (28-30  May,  1984). 


THE  PREDICTION  OF  FORCE  COEFFICIENTS  FOR  LABYRINTH  SEALS 


Otto  W. K.  Lee,  M.  Martinez-Sanchez  and  Eva  Czajkowski 
Massachusetts  Institute  of  Technology 
Cambridge,  MA  02139,  U.S.A. 


1.  Introduction 


Fluid -dynamic  forces  arising  from  nonuniform  pressure  patterns  in  labyrinth  seal 
glands  are  known  to  be  potentially  destabilizing  in  high  power  turbomachinery.  A 
well  documented  case  in  point  is  that  of  the  space  Shuttle  Main  Engine  turbopumps 
(ref.  1),  and  other  examples  can  be  found  in  the  literature,  as  for  instance  in  the 
recent  review  of  Ehrich  (ref.  2)  and  Childs  and  Ehrich  (ref.  3).  Seal  forces  are  al¬ 
so  an  important  factor  for  the  stability  of  shrouded  turbines,  acting  in  that  case  in 
conjunction  with  the  effects  of  blade-tip  clearance  variations  (refs.  4,5). 

The  basic  mechanisms  which  produce  the  uneven  pressure  distribution  in  a  laby¬ 
rinth  have  been  qualitatively  or  semi-qualitatively  discussed  in  many  references 
(refs.  6,7).  In  most  instances,  the  basic  agent  is  found  to  be  flow  swirl  in  the 
glands,  either  from  pre-swirl  (as  in  the  case  of  turbine  shrouds)  or  from  frictional 
interaction  with  the  rotating  shaft  (as  in  multichamber  jet  engine  seals).  Quantita¬ 
tive  modeling  of  these  forces  has  also  been  reported  by  a  number  of  authors  (refs.  8, 
9,10),  using  lumped-parameter  models  for  each  gland.  These  models  yield  in  general 
predictions  of  the  direct  and  cross-wise  stiffnesses  and  damping  coefficients  for 
small  shaft  displacements,  and  are  useful  for  linear  stability  analyses.  Non-linear 
predictions  for  fully  developed  unstable  operation  are  less  advanced. 

Fairly  extensive  data  also  exist  on  the  stiffness  factors  of  seals  of  various  geo¬ 
metries  (refs.  4,5,11,12,13).  These  have  been  generally  obtained  in  rotary  rigs  with 
adjustable  shaft  eccentricity.  Much  less  satisfactory  is  the  situation  with  respect 
to  data  on  damping  coefficients  due  to  labyrinth  seals,  since  these  require  dynamic 
measurements  on  either  vibrating  shafts,  or  shafts  fitted  with  adjustable  whirl  mech¬ 
anisms.  Yet  these  data  are  almost  as  essential  as  those  on  stiffnesses,  since  the 
corresponding  induced  forces  are  of  the  same  order.  Ref. (13)  reports  damping  data 
for  non-rotating  shafts. 

In  this  paper  we  report  on  the  development  of  a  linear  model  for  the  prediction 
of  labyrinth  seal  forces  and  on  its  comparison  to  available  stiffness  data.  We  also 
present  a  discussion  of  the  relevance  of  fluid  damping  forces  and  report  on  the  pre¬ 
liminary  stages  of  a  program  to  obtain  data  on  these  forces. 

2.  Model  Formulation 

The  model  is  very  similar  in  its  main  outline  to  those  of  Kostyuk  (ref.  8)  and 
Iwatsubo  (ref.  9).  It  describes  the  flow  of  an  ideal  gas  through  the  seal  chambers, 
assuming  largely  constant  temperature,  but  allowing  for  isentropic  acceleration  to¬ 
wards  the  narrow  gaps  and  also  for  isentropic  azimuthal  flow  redistribution  in  each 
chamber.  Each  chamber  is  assigned  a  pressure  P.  and  azimuthal  velocity  c^,  and  these 
quantities  are  governed  by  equations  of  mass  ana  azimuthal  momentum  conservation, 
written  in  integral  form. 

The  axial  flow  rate  q  through  each  seal  throttling  is  approximated  by  a  commonly 
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used  expression,  which  basically  derives  from  Bernoulli's  equation  with  a  density 
halfway  between  those  of  the  two  adjacent  chambers.  Per  unit  length  in  the  azimuthal 
direction,  Chis  gives 

P*  -  p2 

i-1  Pi 


where  6^  is  the  width  of  the  narrow  gap  (see  fig.  1)  and  -  c  0^  is  the  product  of 
the  usual  contraction  coefficient  c  times  a  carryover  factor  0  l  to  account  for 
nonzero  upstream  axial  velocity  and  nonzero  pressure  recovery  in  the  downstream  cham¬ 
ber.  Eq.  (1)  is  assumed  to  apply  locally  at  each  time  t  and  azimuth  $.  This  semi- 
incompressible  approximation  is  known  to  be  reasonable  up  to  gap  Mach  numbers  of 
about  0.5;  however,  the  last  gap  or  two  of  a  labyrinth  with  a  high  overall  pressure 
ratio  may  be  above  that  M  ach  number,  and,  in  particular,  the  last  chamber  may  choke. 
We  have  partially  accommodated  this  effect  by  retaining  Eq.  (1)  throughout,  but  re¬ 
placing  it  by  a  choked-flow  expression  in  the  last  chamber  only  if  the  first  approxi¬ 
mation  indicates  sonic  or  supersonic  conditions  there. 

With  reference  to  the  geometry  of  fig.  1,  the  governing  equations  within  each 
chamber  are 


3c  Cplfl3  +  3w  £pificl^  +  qi+l 


3c  (ol£lci)  +  TS  (pl£lc(>  +  "i+l  Wi-1  +  Ti“'  -  T?"  +  fi  TT  *  0  (3) 

These  equations  are  first  linearized  about  a  condition  of  zero  eccentricity.  The 
zero' th  order  approximation  provides  a  basic  flow  rate  q  and  pressure  and  azimuthal 
velocity  distributions  P*  ,  c*  ,  (Appendix  1) .  The  first  approximation  then  provides 
linear  equations  for  the  perturbations,  defined  by 

p.^  ■  p*  (1  +  ;  qi  ■  q*  (1  +  ;  ci  »  ci  (1  +  n1>  (4) 

where  G^»  are  functions  of  t  and  w  »  Rgd>  .  The  right  hand  sides  of  these 
equations  are  determined  by  an  assumed  eccentric  motion  of  the  shaft,  whose  center 
follows  an  elliptic  path 

x  »  r  cos  f2t  ;  y  ■  r  sin  fit  (5) 

Cl  C  2 

where  fi  is  the  shaft  vibration  frequency,  closely  identified  with  one  of  its  natural 
frequencies. 

The  details  of  the  analysis  are  given  in  references  14  and  15.  For  convenient 
solution,  the  perturbations  (for  a  stationary  oscillation)  are  expressed  in  the  form 

£  -  Re  [e  uc(£ssin  $  +  C^cos  $)]  ,  etc.  (6) 


14 


where  R  indicates  the  real  part,  and  £g>  ££  are  in  general  complex  numbers.  After 
elimination  of  5  ,  the  following  system  of  perturbation  equations  is  obtained: 
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The  expressions  for  the  coefficients  A^  ,  B  .  Z^  are  given  in  Appendix  2. 

For  a  seal  with  K  chambers,  i  would  range  from  i  *  1  to  K.  The  calculations  reported 
here  have  assumed  uniform  inlet  and  exit  conditions. 

* 

The  forces  are  then  obtained  by  integration  of  the  perturbation  pressures 
around  the  periphery  of  each  chamber,  followed  by  summation  for  all  chambers.  Since 
our  attention  will  be  restricted  to  the  practical  case  of  circular  whirl  (ri*  r2*  r), 
it  is  advantageous  to  project  the  forces  in  the  directions  towards  the  instantaneous 
minimum  gap  (Fd>  and  90°  ahead  of  it  in  the  whirl  direction  (F^),  as  shown  in  fig.  2. 
We  obtain 

2Tr  * 
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27r  * 

F  -  -R  Z,  f  P.£.  sin  (<fr  -  fit)  d<p  (8b) 
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These  components  can  be  expected  to  be  time- invariant  for  a  symmetric  shaft  in 
circular  whirl.  For  this  case,  the  form  of  the  system  of  equations  (7)  indicates 
that  we  must  have 
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and  the  system  (7)  reduces  to 
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Also,  the  complex  force  vector  F^  +  i  F  is  found  to  be  given  simply  by 


V  + 1  V  *  R»  pi  <- 1 «.  >*t 


(id 


whereas,  if  the  forces  F  and  F  along  fixed  directions  Ox,  Oy  are  desired,  they 
are  given  by  y 

Fx  +  1  Fv  "  (Fd  +  1  Ffl  >  e±nt  (12) 

i  yi  ai  qi 

3.  The  stiffness  and  damping  coefficients 

Most  data  reported  to  date  on  seal  forces  refer  to  situations  with  a  static  off¬ 
set,  i.e.,  with  0*0.  Within  a  linear  approximation,  it  is  then  unambiguous  to  de¬ 
fine  the  direct  and  crosswise  stiffness  K  ,  K  ,  as  -  F,/r  and  -  F  /r  ,  respectively. 
In  practice,  however,  unstable  whirl  is  oisserveS  to  occur  at  or  nea?  (one  of)  the 
shaft  natural  frequencies,  usually  at  the  first  one,  0  (and  it  first  shows  up  at 
rotational  frequencies  co  of  the  order  of  twice  this  natural  frequency).  For  purposes 
of  dynamic  modeling,  then,  it  is  of  interest  to  calculate  or  measure  the  forces  for 
ft  =  ftQ  *  3lnce  t(ie  fluid  disturbances  are  expected  to  translate  into  relatively  small 
real  and  imaginary  departures  of  Q  from  its  basic  resonance  value.  A  common  expres¬ 
sion  for  the  disturbance  forces  in  terms  of  the  x  ,  y  deviations  of  the  shaft  center 
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where  the  stiffnesses  K  and  damping  factors  C  ,  are  taken  to  be  constant.  This 
amounts  to  an  approximation  in  which  terms  proportional  to  x  ,  x  ,  etc.,  are  all 
neglected.  It  is  easy  to  show  that  the  corresponding  approximation  in  terms  of  H  is 


one  where  the  actual  functions  Fd(ft),  F  (“)  are  replaced  by  their  tangents  at  the  fre¬ 
quency  of  interest:  q 
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and  then  the  K  and  C  coefficients  are  given  by 
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Thus,  determination  (either  analytical  or  experimental)  of  F  ,  F  at  two  whirl¬ 
ing  frequencies  near  5)  is  sufficient  to  extract  the  K  and  C  coefricieSts  in  this 
formulation.  ° 

An  alternative  formulation  is  often  found  in  the  literature  (refs.  8,9,15)  and 
is,  in  fact,  the  basis  for  the  results  presented  here.  Analytically,  it  consists  of 
solving  the  system  of  equations  (7)  for  r2  «  0,  i.e.,  for  a  linear  vibratory  shaft 
motion,  and  adscrlbing  to  -K  x  ,  K  x  those  forces  F.  ,  F  in  time  phase  with  x  , 
while  attributing  to  -C  x  ,  -C  x  y  those  in  quadrature  4  with  x  .  It  can  be 
shown  that  this  leads  xx  c  xy  c  tQ  definitions  c 
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A  geometrical  interpretation  of  the  difference  between  equations  (15)  and  (16) 
is  shown  in  fig.  3.  Experimental  or  analytical  determination  of  the  set  (16)  of  co¬ 
efficients  requires  data  on  F.  and  F  at  both  0  and  -fl. 

d  q 

The  K's  and  C's  given  by  (15)  and  (16)  coincide  only  if  F^  and  F  are  linear 
functions  of  (2.  Since  either  definition  may  be  used  in  the  q literature 

with  little  elaboration,  it  would  be  of  interest  to  study  the  extent  to  which  this 
leads  to  numerical  differences;  pending  this,  we  will  in  this  study  adopt  the  defini¬ 
tions  (16).  An  example  for  a  single-chamber  seal  is  shown  in  Appendix  3.  From  this 
limited  evidence,  it  appears  that  the  important  coefficients  K  ,  C  are  about  the 
same  in  both  definitions.  501 


4. 


Comparison  to  Literature  Data 

4. 1  Data  Used 

Benckart  and  Wachter  have  published  an  extensive  set  of  data  (refs.  11,16)  for 
multichamber  labyrinth  seals  of  simple  "straight"  or  "full"  types  (fig.  4).  The 
data  were  taken  in  a  static-offset  rig  operating  on  pressurized  air,  and  induced 
forces  were  obtained  by  integration  of  measured  azimuthal  pressure  variations  on  a 
number  of  seal  cavities.  Labyrinths  with  up  to  23  chambers  were  used.  The  experi¬ 
ments  allowed  variation  of  shaft  speed,  co,  overall  pressure  ratio  P  /P  ,  rotor 
eccentricity  r,  seal  geometry  (6i  ,  ,  h  ,  Fig.  1),  number  of  atm  0  chambers 

and  entry  swirl  c*  .  The  seal  flow  rate  q  was  measured  and  an  averaged  carry¬ 
over  factor  u  wa S  deduced  from  these  data  and  reported  in  a  number  of  instances.  In 
our  calculations  we  used  these  "measured"  factors  when  available  directly;  in  other 
cases,  we  adopted  values  measured  for  chambers  of  the  same  geometry,  or,  for  the 
"full"  type  of  seal,  where  little  carryover  is  expected  we  used  u  -  1.  The  contrac¬ 
tion  coefficient  c  was  taken  as  a  function  of  Reynolds  number  and  strip  geometry  as 
given  by  Vermes  (ref.  17)  (fig.  5). 

Data  of  Brown  and  Leong  (ref.  18)  were  also  used  for  validation  of  the  undis¬ 
turbed  flow  predictions. 

4.2  Undisturbed  Flow  Parameters 


Figure  6  compares  Brown  and  Leong 's  data  on  the  axial  pressure  distribution  in 
an  11-chamber  test  seal  with  our  calculated  undisturbed  pressure  distribution.  There 
is  good  agreement  except  for  the  sharp  pressure  drop  shown  by  the  data  between  the 
inlet  and  the  first  chamber.  This  is  probably  a  reflection  of  a  reduced  carry-over 
factor  on  the  first  strip;  the  calculation  used  a  constant  u  (the  value  used  is  ir¬ 
relevant  to  the  comparison) . 

Benckert  and  Wachter  reported  for  one  particular  case  the  axial  variation  of 
azimuthal  velocities  c*  .  This  was  for  a  23-chamber  seal  with  c*  ■  40  m/sec, 

Pa/PQ  ■  .66,  Rg  *  0.15  m,  <jj  -  1000  rad/sec,  n  ■  0.25  mm,  T2  m  0°  6*  *  0.5  mm, 

H  *  4  mm  and  h^  ■  6  mm.  The  data  are  shown  in  figure  7,  together  with  the  code  pre¬ 
dictions.  The  good  agreement  shown  is  important  for  the  prediction  of  disturbance 
side  forces,  which  depend  critically  on  swirl  velocities.  These  results  appear  to 
validate  the  formulation  used  for  the  friction  factors  between  the  fluid  and  the 
stator  and  rotor  surfaces  (turbulent  pipe  flow  formulae  with  a  standard  correction 
for  "pipe"  curvature). 

4.3  Stiffness  Coefficients  Without  Shaft  Rotation 


The  cross-spring  coefficients  K  -  -  K  for  a  number  of  cases  from  Benckert 
and  Wachter' s  tests  with  a  non-rotating  shaft^were  calculated  and  the  results  are 
summarized  in  Table  2  and  figure  8.  The  key  in  Table  2  describing  the  test  para¬ 
meters  is  explained  in  Table  1. 

The  eccentricity  ri  used  in  the  tests  was  0.15  mm,  except  for  Run  17,  which  had 
n  *  0.25  mm. 


Table  1.  Key  for  Table  2  (2nd  column) 


S  ■  Straight- through 
F  »  Full- interlocking 


(Type,  h,  6,  Z) 


h:  chamber  height 


0  ■  2.75  mm 
1  -  6.25  mm 


clearance 


■  0.25  mm 
■0.5  mm 


Z:  pitch 


As  figure  8  illustrates,  the  calculated  values  are  somewhat  lower  than  the  data 
(about  19%  for  series  (s, 0,0,0),  5%  only  for  series  (s, 1,0,0)).  The  trends  of  the 
calculation  are  in  agreement  with  those  observed  in  the  tests.  In  particular,  K 
is  seen  in  figure  9  to  be  approximately  proportional  to  inlet  swirl  and  to  overafi 
pressure  ratio,  both  in  the  tests  and  in  the  calculations  (although,  as  mentioned, 
with  a  somewhat  lower  proportionality  factor  in  the  latter  case). 

4.4  Stiffness  Coefficients  with  Shaft  Rotation 


Results  similar  to  those  in  the  previous  section,  but  including  various  shaft 
speeds  are  given  in  Table  3  and  figure  10.  The  eccentricity  is  still  static  (no 
whirl,  0  ■  0) ,  and  is  0.25  mm  in  all  cases.  The  parameter  E*  was  used  in  Ref.  11  to 
correlate  entry  swirl,  and  is 

.  (c*)2  * 

E  -  - f-2 - - -  ;  V  -  (17) 

O  .  1  „2  ax,0  x* 

p  V  -  p  p  Si 
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The  comparison  of  data  and  theory  shown  in  figure  10  indicates  more  scatter,  but 
less  systematic  deviation  than  in  the  cases  without  shaft  rotation  (figures  8  and  9). 
The  agreement  is  best  for  all  the  cases  with  17  chambers  (solid  symbols  in  figure  10), 
which  show  an  average  error  of  8.5%  and  little  scatter. 

4.5  Discussion 


The  two  principal  sources  of  uncertainty  in  our  calculations  are  the  friction 
factors  (X',  X")  and  the  carry-over  coefficient  S.  The  friction  factor  could  in 
principle  be  substantially  increased  by  the  relative  rotation  of  shaft  and  casing, 
since  the  fluid  in  each  chamber  is  strongly  sheared  and  develops  marked  secondary 
flow  patterns,  leading  to  enhanced  mixing.  Examination  of  data13for  the  somewhat  re¬ 
lated  case  of  turbulent  pipe  flow  with  swirl  does  indeed  show  friction  increases  of 
up  to  a  factor  of  four  at  high  swirl.  An  accurate  prediction  of  wall  friction  under 
the  complex  flow  conditions  of  a  labyrinth  gland  is  not  possible  at  this  time,  and 
this  is  an  area  requiring  more  experimental  and  analytical  work. 

The  impact  of  friction  factor  inaccuracies  on  calculated  cross-spring  coeffici¬ 
ents  could  be  important,  although  not  easily  generalizable.  In  general,  the  cross¬ 
forces  increase  with  the  deviation  between  the  swirl  velocity  and  its  frictional 
equilibrium  value.  An  increase  in  wall  friction  in  a  non-rotating  seal  will  acceler¬ 
ate  the  approach  to  this  ultimate  swirl,  thus  reducing  the  number  of  chambers  where 
the  excess  or  defect  swirl  is  strong  and  thereby  reducing  the  magnitude  of  Kyx 


(whether  positive  or  negative).  This  trend  is  always  apparent  at  large  enough  fric¬ 
tion  coefficients,  and  especially  for  long  multichamber  seals.  However,  for  short 
seals  and  weak  frictional  coupling  a  different  effect  dominates,  namely  the  cross¬ 
stiffness  tends  to  a  limiting  value  independent  of  shaft  rotation  and  in  the  direc¬ 
tion  of  the  inlet  swirl.  An  example  of  this  behavior  for  a  one-chamber  seal  is  pre¬ 
sented  in  Appendix  3,  figure  A2,  where  it  can  be  seen  that  for  this  particular  case, 
increasing  friction  would  lead  to  increases  in  |k  |.  The  behavior  typical  of  long 
seals  and  large  friction  is  illustrated  in  figure^*  11,  corresponding  to  case  1  of 
Table  2;  here,  an  increase  of  the  friction  factor  above  the  nominal  value  leads  to  a 
reduced  K  ,  although  a  reduction  by  more  than  about  0.6  would  lead  to  the  same 
effect. 

The  carry-over  coefficient  8  is  clearly  another  source  of  uncertainty  in  the 
model.  The  sensitivity  of  calculated  cross-stiffness  to  8  (or  u  ■  8c  )  is  illustra¬ 
ted  in  figure  12,  corresponding  to  parametric  variations  on  Run  11  ofcTable  2,  and 
figure  13  (from  Run  7  of  Table  3).  The  opposite  trend  in  these  two  cases  is  due  to 
the  fact  that  in  figure  12  the  entry  swirl  is  greater  than  its  asymptotic  value 
(reached  after  an  infinite  number  of  chambers),  while  the  opposite  is  true  in  figure 
13.  In  both  cases,  an  Increase  of  u  increases  the  flow  rate  q,  which  has  the  effect 
of  delaying  the  transition  towards  the  asymptotic  c*  ;  in  figure  12  this  means  higher 
c*  in  the  first  10-12  chambers,  with  correspondingly . larger  cross-forces;  in  figure 
13,  the  same  delaying  effect  at  higher  u  implies  lower  c*  values  in  the  first  6-8 
chambers,  and,  consequently,  lower  cross-forces. 

Unfortunately,  the  state  of  the  art  in  a  priori  predictions  of  8  is  not  satis¬ 
factory.  A  fuller  discussion  of  this  point  is  given  in  reference  15.  Basically,  the 
best-known  models  (Vermes  (ref.  17),  Egli  (ref.  20),  Komotori  (ref.  21)),  indicate 
substantially  different  variations  of  the  carry-over  factor  with  number  of  chambers, 
with  Vermes'  model  taking  no  account  of  this  number  at  all.  The  best  hope  here  lies 
with  the  numerical  methods  which  are  now  beginning  to  be  applied  to  internal  flow 
problems  in  seals,  although  the  somewhat  primitive  state  of  affairs  with  respect  to 
calculations  of  fully  separated  turbulent  flows  still  indicates  a  need  for  improve¬ 
ments.  Thus  Wadia  and  Booth  (ref.  22)  analyzed  seal  flows  with  no  rotation  and  ob¬ 
served  discrepancies  of  up  to  13%  in  calculated  flow  coefficients  when  compared  to 
data.  For  dynamics  studies  in  seals,  these  2-D  or  3-D  methods  may,  in  any  case,  be 
too  laborious;  their  proper  role  should  probably  be  in  furnishing  improved  semi- 
empirical  results  for  integration  into  a  simple  multi-chamber  lumped-parameter  model, 
of  the  type  considered  here. 

5.  Parametric  Studies 

Reference  15  includes  a  variety  of  calculations  that  illustrate  the  trends  of 
the  force  coefficients  versus  variations  of  seal  parameters.  Only  some  of  the  sal¬ 
ient  results  will  be  mentioned  here. 


(a)  K  increases  linearly  (but  not  proportionally)  with  entry  swirl  velocity.  For 

conditions  where  the  entry  swirl  exceeds  the  asymptotic  azimuthal  velocity,  K  is 

generally  positive,  leading  to  excitation  of  forward  whirl  (with  respect  to  the 

swirl  direction).  The  reverse  may  be  true  at  lower  entry  swirls.  There  is  in  some 

cases  a  value  of  entry  swirl  at  which  K  is  zero. 

yx 

(b)  For  multichamber  seals  with  low  entry  swirl,  the  first  few  chambers  contribute 

negative  K  values,  while  those  towards  the  end  of  the  seal  contribute  positive 

values.  Thus  as  the  seal  is  made  longer,  the  sign  of  K  may  at  some  point 


reverse.  For  seals  where  the  entry  swirl  exceeds  the  ultimate  azimuthal  velocity,  no 
such  reversal  occurs. 


(c)  The  damping  coefficient  C  ,  which,  together  with  K  controls  the  side  force  F  , 
is  positive  in  all  cases  studied,  leading  to  stabilizing^^orces  of  the  same  order  as^ 
those  due  to  K  .  This  point  will  be  more  fully  discussed  in  what  follows. 

y* 

6.  Considerations  on  Fluid  Damping 

An  example  of  calculated  damping  coefficients  is  shown  in  figure  14  (correspond¬ 
ing  to  the  seal  configuration  of  Run  1  in  Table  3).  A  whirl  (critical)  frequency  of 
739  rad/sec  was  assumed;  at  the  commonly  found  ratio  w/ft  »  2  for  instability  onset, 

C  ■  220  N  sec/m  ,  giving  SlCxx  ■  1. 626x10 sN/m.  This  is  several  times  larger  than 

K**,  and  indicates  that  seal  damping  forces  are  in  this  case  sufficient  to  en- 

^  sure  stability.  Another  example  of  this  behavior  is  shown  by  the  single-cavity 
seal  of  Appendix  3;  here  +  is  roughly  comparable  to  K^x  (but  smaller). 

A  very  simplified  dynamic  model  will  help  to  put  in  perspective  the  roles  of  the 
different  coefficients  in  stability  analysis.  Assuming  a  shaft  with  mass  M  and  struc¬ 
tural  stiffness  Kq,  the  equations  of  motion  for  small  side  displacements  x,y  can  be 
combined  into 

A  A 

Mz+Cz+Kz-0  (18) 

A  A 

where  C  ■  C  +  1  C  ,K»K+K  +  i  K  and  z  *  x  +  iy. 
xx  xy  o  xx  xy 

Assuming  K  ,  K  ,  C2  /4M  and  C2  /4M  are  all  small  compared  to  K  ,  as  is  likely  to  be 
8  xx  xy  xx  xy  ° 

the  case  in  practice,  we  can  define  the  (small)  nondimens ional  parameters 


2i/K~M 

o 


(19) 


Then  a  simple  analysis  shows  that, 


to  the  first  approximation,  the  shaft  complex 


displacement  3  will  vary  as  e 


figt 


where 


A  /M 
o 


£  iU  )+i(±l+^k  -  £  ). 

xx  2  xy  2  xx  xy 


(20) 


Thus  any  nonzero  k  will  be  destabilizing  (in  one  or  the  other  whirl  direction), 

while  a  negative  £  wi¥i  be  always  destabilizing;  k  and  £  will  simply  modify 

the  shaft  natural  ^frequency.  Also,  the  effects  of  equal  va^taes  of  |K  |  and 

-  ft  C  are  seen  to  be  equivalent.  “y 

xx 


This  discussion  has  served  to  indicate  that  knowledge  of  the  damping  factor  C xx 
is  at  least  as  essential  to  studies  of  fluid-induced  destabilizing  forces  as  is 
the  side  force  factor  K  .  Yet,  due  to  the  more  difficult  experimental  conditions, 
much  fewer  data  are  avaHable  on  C  than  on  K 


We  are  now  in  the  design  stage  of  a  test  rig  intended  to  address  this  problem. 

The  general  size  and  flow  parameters  will  be  similar  to  those  used  by  Benckert  and 
Wachter,  but  the  sealed  shaft  will  be  made  to  execute  forced  whirling  motion  at 
speeds  controlled  separately  from  the  spinning  motion.  Pressure  distributions  will 
be  dynamically  measured  and  integrated  to  produce  values  of  the  direct  and  transverse 
forces  F,  ,  F  for  a  range  of  whirl  speeds  £).  Both  K  and  C  coefficients  can  then  be 
extracted  by  ?he  methods  described  in  this  paper.  These  features  are  similar  to  those 
of  a  rig  described  in  reference  23  for  tests  in  water. 


7.  Conclusions 


A  linear  analytical  model  for  the  prediction  of  fluid  forces  in  labyrinth  seals 
has  been  presented  and  discussed.  Comparison  to  literature  test  data  shows  reason¬ 
able  agreement  for  the  important  cross-stiffness  K  .  The  importance  of  the  damping 
factor  C  has  been  highlighted  and  the  need  for  damping  data  made  clear. 


xx 
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Appendix  1.  The  Unperturbed  Solution 

Squaring  Eq.  (1)  and  adding  for  all  chambers  yields  for  the  nominal  flow  rate 


P2  -  p2 

o  a 


(A-l) 


RTS  (- 


“i  s‘i! 


m 


J:::;;:' 

•e  V  ..J 


W2 

•v-v-v 


-\V 

vC-KO' 

fr  4 


•  V  /  4 

*v  V  „ 
\  v  «. 

"  V.  , 


Also,  adding  for  the  first  n  chambers  only  gives 


A  ,  A 

p  2  ,  p  2 

n  o 


Ji  W 

K 


(P2  -  P2) 
o  a 


Z  <-T 

i-1  uJV 

The  momentum  equation  (Eq.  (3))  becomes  in  the  steady  state 


*  * 


q  (c1  -  ci_1)  +  Tj  U'  -  Tj  U"  -  0 


(A-2) 


(A- 3) 


with 


Ti  ■  i p* x'  4* 


(A-4) 


T"  -  i  P*  X"  <0*  -  C*)! 


(A-5) 


and  X  (the  Darcy  friction  factor)  given  by  a  modified  pipe-flow  expression 


,  0.3164  r,  .  „  „„  „  0.25,Dh,1/2,  .  ,, 

X  *  ~OJ  11  +  °-075  Re  (2r}  1  x  3ign  (VREL)  (A_6) 

t\  3 

e 

Here  R  is  the  Reynolds  number  based  on  chamber  height  and  the  corresponding  relative 
flow  velocity  for  fixed  or  rotating  surfaces.  The  factor  sign  (v  )  is  needed 
to  give  the  forces  their  proper  direction.  Thus,  we  append  tne  factor  sign 

(c*)  to  X'  and  the  factor  sign  (oiR  -  c*)  to  X"  ,  both  here  and  in  the  first  order 
calculations  of  Appendix  2.  is3  the  gland  hydraulic  diameter.  Eqs.  (A-3) 
through  (A-6)  can  be  solved  for  the  distribution  c.  of  azimuthal  velocities.  In 
particular,  the  asymptotic  velocity  (cot  )  follows  from  (A-3)  when  c*  *  c*_i  * s 
assumed. 


Appendix  2.  Coefficients  for  the  Perturbation  Equations 


★  *  * 

Let  the  transverse  area  of  a  gland  (fig.  1)  be  f  *  (h  +  6  )  in  the  centered 
position.  An  asterisk  on  any  variable  denotes  the  undisturbed  (centered)  condition. 
We  obtain  for  Eq.  (7)  the  following  coefficients. 


A 


i 


★  ★ 


Rs  Y 


* 

Ci  f 


s 


* 

i 


n  - 


Pi+1  Ui  5i+l 


F. 


* 

P 


U 


2 

i 


i+1  i+1  i+1 


i  ^i+1  i+1 


U*  X'  c 


-izi 

*  * 


W  -  a_  (  _ izi  -J- 

"j  *  v  *  *  * 


*i  ■  *i 


* 

ci  h 


Here  y  is  the  ratio  of  specific  heats  and  the  gas  constant. 


Appendix  3.  "Local"  vs.  "Global"  Coefficients 


Sample  calculations  were  made  for  a  single-chamber  straight-through  seal  with 
6  *  0.25  mm.,  2,  ■  8  mm.,  h  ■  3  mm.,  R  *  15  cm.,  w  ■  1000  rad/sec,  c  ”  100  m/sec , 

P  ■  1.5  atm,  P  *  1  atm  .  For  a  ran|e  of  whirl  frequencies  from  ft  »  0  to  ft  *  750 
0  rad/sec,  the3  resulting  direct  and  quadrature  forces  are  shown  in  figure  Al.  The 
quadrature  force  F  ,  which  is  the  one  of  importance  for  stability  considerations,  is 
seen  to  be  very  ^  nearly  linear  with  ft,  indicating  no  difference  between  the  local 
and  the  global  definitions  of  the  coefficients  (Eqs.  (15)  and  (16),  respectively). 


There  is, 
lated  for 

on  the  other  hand,  a  slight  curvature 
ft  *  500  rad/sec  are  as  follows: 

in  the  F,  line, 
a 

The  values  calcu 

K  (N/m) 

XX 

K  (N/m) 
yx 

C  (N  sec/m) 

XX 

C  (N  sec/m) 
yx 

Local 

(Eqs.  (15)) 

610 

11350 

+  13.7 

-  0.95 

Global 

(Eqs .  (16)) 

652 

11480 

13.98 

-  1,02 

For  the  same  seal,  with  ft  =  0  throughout,  figure  A2  shows  the  effect  of  para¬ 
metric  variations  of  the  friction  coefficients  (V  and  X"  varied  simultaneously)  at 
various  shaft  rotation  speeds. 


K  Calculated  vs.  Experimental  (uj^0 t 


Run  #  Seal  Type 


(x  10  N/m) 


l 

S, 0,0,0 

17 

.66 

38.4 

.92 

.75 

.611 

2 

S,  0,0,0 

17 

.32 

52.9 

.92 

2.57 

2.091 

3 

S, 1,0,0 

17 

.32 

68.3 

1.02 

1.57 

1.45 

4 

S, 1,0,0 

17 

.66 

33.4 

1.02 

.27 

0.274 

5 

S, 1,0,0 

17 

.66 

48.5 

1.02 

.423 

0.411 

6  . 

S, 1,0,0 

17 

.56 

39.0 

1.02 

.457 

0.410 

7 

S, 1,0,0 

17 

.  79 

38.2 

1.02 

.218 

0.224 

8 

S, 0,0,0 

17 

.49 

63.5 

.92 

1.89 

1.586 

9 

S, 0,0,0 

17 

.39 

54.3 

.92 

2.22 

1.747 

10 

S, 0,0,0 

17 

.79 

15.6 

.92 

.184 

0.160 

11 

S, 0,0,0 

17 

.49 

64.2 

.92 

1.75 

1.605 

12 

S, 0,0,0 

17 

.49 

34.5 

.92 

1.05 

0.807 

13 

S, 1,0,0 

17 

.49 

82.68 

1.02 

.98 

1.098 

14 

S, 1,0.0 

17 

.49 

40.78 

1.02 

.57 

0.509 

15 

S, 0,0,0 

17 

.  32 

38.2 

.92 

1.9 

1.431 

16 

S, 0.0,0 

17 

.  32 

27.6 

.92 

1.2 

0.937 

17 

F, 3,2,1 

9 

.49 

144.7 

.665 

1.47 

1.606 

18 

S, 0,0,0 

17 

.  32 

45.5 

.92 

2.39 

1.757 

Average  error  (in 

absolute 

value)  * 

18. 32 

(4.52 

for  (S,l,0, 

0)  , 

,  18.62 

for  (S.O 

.0,0)) . 

Table 

3.  K  Calculated  vs. 
y* 

Experimental  (wrot 

*  0) 

K 

y* 

P  f? 
a  o 

0 

rot 

Exper.  gCalcula 
(x  10  N/m) 

Run  # 

Seal  Type 

K 

(a/s) 

(rad/s) 

u 

F.3,2,2 
F.3.2,2 
F.3,2,2 
F.3.2.2 
F. 3,2,2 
F.3,2,2 
F.3,2,2 
F.3,2,2 
F.3,2,2 
F. 3,2,2 
F.3,2,2 


746.67 

993.3 

993.3 

746.67 

1000 

1000 

746.67 

500. 

1000 

746.67 

500 


0.228 
0.230 
0.  398 
0.  340 
0.088 
0. 187 
0.0849 
0.0178 
0. 328 
0.  170 
0.123 


1  —  U  calculated  from  measured  mass  flow  rate:  friction  factor  for  a 

channel  used. 

2  —  resultant  forces  from  chambers  7-17  only. 

3  --  resultant  forces  from  chambers  7-23  only. 

Average  error  (in  absolute  value)  *  23.02 


Fig.  1  Geometry  for  labyrinth  seal  analysis 


Fig.  2  Direct  (F  )  and  quadrature  (F  )  forces 
P  <? 

due  to  rotor  eccentricity 


Fig.  3  Two  different  definitions  of  spring  and 
damping  coefficients. 

(a)  Local  (Eqs.15),  (b)  Global  (Eqs.16)) 


Coefficient  of  contraction,  C 


Ky x  CAL 


Kyx  CALCULATION 


M 


Kyx  EXPERIMENT  (  X  t0~5  N/M) 

Fig.  10.  Calculated  va.  experimental  croaa  spring  coefficient  Kyx  (A>rot  i  0) 


FRICTION  FACTOR  RATIO  V^cp 


Fig.  11.  Effect  of  friction  factor  on  cross  spring  coefficient  Kyx. 
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NON  DIMENSIONAL  ROTATION  FREQ. 

fflclents  va.  non-dlaenaional  rotation  frequency 


